﻿using MathNet.Numerics;
using System;
using System.Collections.Generic;
using System.Linq;

namespace Jodell.Lua
{
    public class ParticleSwarmOptimization
    {
        /// <summary>
        /// 粒子群大小
        /// </summary>
        private int swarmSize;
        /// <summary>
        /// 解空间的维度
        /// </summary>
        private int dimensions;
        /// <summary>
        /// 惯性权重
        /// </summary>
        private double inertiaWeight;
        /// <summary>
        /// 认知权重
        /// </summary>
        private double cognitiveWeight;
        /// <summary>
        /// 社会权重
        /// </summary>
        private double socialWeight;
        /// <summary>
        /// 解空间的下界
        /// </summary>
        private double[] lowerBounds;
        /// <summary>
        /// 解空间的上界
        /// </summary>
        private double[] upperBounds;
        /// <summary>
        /// 全局最优适应度
        /// </summary>
        public double GlobalBestFitness { get; private set; }
        /// <summary>
        /// 全局最优位置
        /// </summary>
        public double[] GlobalBestPosition { get; private set; }
        /// <summary>
        /// 粒子群
        /// </summary>
        private Particle[] particles;
        /// <summary>
        /// 
        /// </summary>
        /// <param name="swarmSize">粒子群大小</param>
        /// <param name="dimensions">解空间的维度</param>
        /// <param name="inertiaWeight">惯性权重</param>
        /// <param name="cognitiveWeight">认知权重</param>
        /// <param name="socialWeight">社会权重</param>
        /// <param name="lowerBounds">解空间的下界</param>
        /// <param name="upperBounds">解空间的上界</param>
        public ParticleSwarmOptimization(int swarmSize, int dimensions, double inertiaWeight, double cognitiveWeight, double socialWeight, double[] lowerBounds, double[] upperBounds, double[] xds, double[] yds)
        {
            this.swarmSize = swarmSize;
            this.dimensions = dimensions;
            this.inertiaWeight = inertiaWeight;
            this.cognitiveWeight = cognitiveWeight;
            this.socialWeight = socialWeight;
            this.lowerBounds = lowerBounds;
            this.upperBounds = upperBounds;
            this.GlobalBestFitness = double.MaxValue;
            this.GlobalBestPosition = new double[dimensions];

            this.XDatas = xds;
            this.YDatas = yds;

            particles = new Particle[swarmSize];

            for (int i = 0; i < swarmSize; i++)
            {
                particles[i] = new Particle(dimensions);
                for (int j = 0; j < dimensions; j++)
                {
                    particles[i].Position[j] = RandomDouble(lowerBounds[j], upperBounds[j]); // 随机初始化粒子位置
                    particles[i].Velocity[j] = RandomDouble(lowerBounds[j], upperBounds[j]); // 随机初始化粒子速度
                    particles[i].BestPosition[j] = particles[i].Position[j]; // 初始化粒子的历史最优位置为当前位置
                }
                particles[i].BestFitness = ObjectiveFunction(particles[i].Position); // 计算粒子的适应度
                if (particles[i].BestFitness < GlobalBestFitness)
                {
                    GlobalBestFitness = particles[i].BestFitness;
                    particles[i].BestPosition.CopyTo(GlobalBestPosition, 0); // 更新全局最优位置
                    globalParticle = particles[i];
                }
            }
        }


        public double[] YDatas;
        public double[] XDatas;
        public Particle globalParticle;
        public void Optimize(int maxIterations)
        {
            for (int iteration = 0; iteration < maxIterations; iteration++)
            {
                //Console.WriteLine("第{0}次迭代",iteration);

                for (int i = 0; i < swarmSize; i++)
                {
                    for (int j = 0; j < dimensions; j++)
                    {
                        double r1 = RandomDouble(0, 1);
                        double r2 = RandomDouble(0, 1);

                        particles[i].Velocity[j] = inertiaWeight * particles[i].Velocity[j] +
                            cognitiveWeight * r1 * (particles[i].BestPosition[j] - particles[i].Position[j]) +
                            socialWeight * r2 * (GlobalBestPosition[j] - particles[i].Position[j]); // 更新粒子速度

                        particles[i].Position[j] += particles[i].Velocity[j]; // 更新粒子位置

                        // 边界处理
                        if (particles[i].Position[j] < lowerBounds[j])
                            particles[i].Position[j] = lowerBounds[j];
                        else if (particles[i].Position[j] > upperBounds[j])
                            particles[i].Position[j] = upperBounds[j];
                    }

                    double fitness = ObjectiveFunction(particles[i].Position); // 计算粒子的适应度
                    if (fitness < particles[i].BestFitness)
                    {
                        particles[i].BestFitness = fitness;
                        particles[i].Position.CopyTo(particles[i].BestPosition, 0); // 更新粒子的历史最优位置
                        if (particles[i].BestFitness < GlobalBestFitness)
                        {
                            GlobalBestFitness = particles[i].BestFitness;
                            particles[i].BestPosition.CopyTo(GlobalBestPosition, 0); // 更新全局最优位置
                            globalParticle = particles[i];
                        }
                    }
                }
            }
        }

        private double RandomDouble(double minValue, double maxValue)
        {
            Random random = new Random();
            return random.NextDouble() * (maxValue - minValue) + minValue; // 生成指定范围内的随机数
        }

        public double ObjectiveFunction(double[] position)
        {

            var xdatas = position.Distinct().ToList();
            xdatas.Sort();
            var min = XDatas.First();
            List<Tuple<double, double>>[] datas = new List<Tuple<double, double>>[xdatas.Count + 1];
            for (int i = 0; i < xdatas.Count; i++)
            {
                for (int j = 0; j < XDatas.Length; j++)
                {
                    if (XDatas[j] >= min)
                    {
                        if (XDatas[j] < xdatas[i])
                        {
                            if (datas[i] == null)
                            {
                                datas[i] = new List<Tuple<double, double>>();
                            }
                            datas[i].Add(new Tuple<double, double>(XDatas[j], YDatas[j]));
                        }
                        else
                        {
                            break;
                        }
                    }
                }
                min = xdatas[i];
            }
            //最后一段
            datas[xdatas.Count] = new List<Tuple<double, double>>();
            for (int j = 0; j < XDatas.Length; j++)
            {
                if (XDatas[j] >= min)
                {
                    datas[xdatas.Count].Add(new Tuple<double, double>(XDatas[j], YDatas[j]));
                }
            }
            int cnt = 0;
            double fitness = 0;
            foreach (var d in datas)
            {
                if (d?.Count > 2)
                {
                    var xs = d.Select(x => x.Item1).ToArray();
                    var ys = d.Select(x => x.Item2).ToArray();
                    var regression = Fit.Line(xs, ys);

                    for (int i = 0; i < d.Count; i++)
                    {
                        var y = d[i].Item2 - (regression.B * d[i].Item1 + regression.A);
                        y = y * y;
                        fitness += y;
                        cnt++;
                    }
                }
            }
            return fitness / cnt;
        }
    }
}
